Young's modulus of Graphene: a molecular dynamics study 
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The Young's modulus of graphene is investigated through the intrinsic thermal vibration in 
graphene which is 'observed' by molecular dynamics, and the results agree very well with the recent 
experiment [Science 321, 385 (2008)]. This method is further applied to show that the Young's 
modulus of graphene: (1). increases with increasing size and saturates after a threshold value of the 
size; (2). increases from 0.95 TPa to 1.1 TPa as temperature increases in the region [100, 500]K; 
(3). is insensitive to the isotopic disorder in the low disorder region (< 5%), and decreases gradually 
after further increasing the disorder percentage. 
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The single layer graphene has unique electronic and 
other physical properties, thus becoming a promising 
candidate for various device applications.^'^ Among oth- 
ers, excellent mechanical property is an important advan- 
tage for the practical appUcations of graphene. Experi- 
mentally, the Young's modulus (Y) of graphene has been 
measured by using atomic force microscope (AFM) to in- 
troduce external strain on graphene and record the force- 
displacement relationi^ The measured value for Young's 
modulus is 1.0 ± 0.1 TPa in this experiment. Theoret- 
ically, the Young's modulus of graphene can be stud- 
ied in a parallel way. Once the external strain is ap- 
plied on graphene, the internal force or potential can 
be calculated in different approaches, such as ab initio 
calculationSf^i^i^ molecular dynamics (MD)Z and inter- 
atomic potentials4i2iiiS Then the Young's modulus can 
be obtained from the force-displacement or the potential- 
displacement relation. For the carbon nanotubes (CNT), 
the Young's modulus is theoretically studied in a similar 
way as that in graphene. However, in the experiment, 
besides the AFM method, another group measured the 
Young's modulus of CNT by observing the thermal vi- 
bration at the tip of the CNT using the transmission elec- 
tron microscopy (TEM) J^iH For some unknown reasons, 
possibly technical challenges, this experimental method 
does not appear in the study of the Young's modulus in 
graphene. As a supplement to this vacancy, the present 
work 'observes' the thermal vibration of graphene by MD 
instead of TEM, and then calculates the Young's modu- 
lus from the 'observed' thermal vibration. 

In the engineering application of graphene, it will be 
beneficial if the mechanical property of graphene can be 
adjusted according to the demand. There are some pos- 
sible methods that can manipulate the value of Young's 
modulus in graphene, such as size of the sample, temper- 
ature, isotopic disorder, etc. It is a matter of practical 
importance and theoretical interest to find an effective 
method to control the mechanical property of graphene. 
The present calculation method for the Young's modulus 
of graphene in this paper is readily applicable to address 
these issues. 



In this paper, we investigate the Young's modulus of 
graphene by 'observing' the thermal vibrations with MD. 
The calculated Young's modulus is in good agreement 
with the recent experimental one. Using this method, we 
can systematically study different effects on the Young's 
modulus: size, temperature and isotopic disorder. It 
shows that the Young's modulus increases as graphene 
size increases, and saturates. In the temperature range 
100 - 500 K, Y increases from 0.95 TPa to 1.1 TPa as 
T increases. For the isotopic disorder effect, Y keeps al- 
most unchanged within low disorder percentage (< 5%), 
and decreases gradually after further increasing disorder 
percentage. 

In graphene there are both optical and acoustic vibra- 
tion modes in the z direction. For the optical phonon 
modes, the frequency is about 850 cm~^, which is too 
high to be considerably excited under 500 K. While the 
acoustic phonon mode is a flexure mode with parabolic 
dispersion ui = (3k'^, which will be fully excited even at 
very low temperature. So the thermal mean-square vi- 
bration amplitude (TMSVA) of graphene in the z direc- 
tion is mainly attributed to the flexure mode under 500 
K. In this sense, we consider the contribution of the flex- 
ure mode to TMSVA for an elastic plate in the following. 
The X and y axes lie in the plate, and z direction is per- 
pendicular to the plate. For convenience and without 
losing generality, we consider a square plate with length 
L. 

The equation for oscillations in z direction of a plate 
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where D — j^Yh^/{l — //^). A is the two-dimensional 
Laplacian and p is the density of the plate. Y and fj, are 
the Young's modulus and the Poisson ratio, respectively. 
h is the thickness of the plate. We apply fixed boundary 
condition in x direction, and periodic boundary condition 
in y direction: 
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z{t,x = L,y) = 0, 
z{t,x,y + L) = z{t,x,y). 



(2) 



The solution for the above partial differential equation 
under these boundary conditions can be found in Ref. il5i : 



Zn{t,x,y) = u„ sin(fcix) • cos(fc22/) • cos(w„t), (3) 
k = kiCx + k2ey, 

where ki — irni/L and /c2 — 27r?i2/L. 

Using these eigen solution, the TMSVA for rt-th 
phonon mode in (x, y) at temperature T can be 
obtained^'^ : 

al{x,y) = 4fcsrx X ■^(sin(fcia;)cos(fc2y){)l: 

We mention that for those modes with fci ^ and fc2 = 0, 
we have a similar result a^{x,y) = 2kBT x ^^j^^f^y ^ x 
^(sin kix)'^. 

The spatial average of the TMSVA over x and y is: 



alix,y)dxdy 
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where ki ^ and k2 0. D is the field in x G [0, i] and 
y G [0, L], and 5* = is the area of D. If /ci 7^ and 
^2 = 0, (cr^) turns out to have the same expression as 
this general one. 

Because all modes are independent at the thermal equi- 
librium state at temperature T, they contribute to the 
TMSVA incoherently. As a result, the TMSVA at tem- 
perature T is given by: 
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7 >^ 1.2507, 



The constant C = ^"^ (;^4n^) 

the major part of which is due to the first nonzero phonon 
mode with (ni,n2) = (1,0). e^^ — 1 for 712 — 0, and 
e„2 — 2 for other 712 = 1, 2, 3, .... 

As a result, the Young's modulus of the graphene is: 
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The Poisson ratio in graphenei^iii fj, = 0.17 has been 
used in this expression for the Young's modulus. There 
is arbitrariness in the definition of thickness h of the one 
atom thick graphene sheet. For convenience of compari- 
son between our theoretical results and the experimental 
ones, we choose h to be 3.35 A, the inter-layer space in 
graphite, which is also used in the experimental worki^ 




FIG. 1: (Color online) Configuration of the graphene sample. 
The origin O is at the left bottom of the sample. Two columns 
(red online) on the left and right sides are fixed. The length 
of the sample in this figure is L = 40 A. 

Fig. [1] is the configuration of the graphene sheet in our 
simulation. The outmost two columns (red online) on the 
left and right sides are fixed during the simulation, while 
periodic boundary condition is imposed in the vertical 
direction. The origin of the coordinate is set at the left 
bottom of the sample. The x-axis is in the horizontal 
direction and y-axis is in the vertical direction. 

The MD simulations are performed using the second- 
generation Brenner inter-atomic potentialJ^ The Newton 
equations of motion are integrated within the fourth or- 
der Runge-Kutta algorithm, in which a time step of 0.5 fs 
is applied. The typical MD simulation steps in this paper 
is 5 X lO"", corresponding to 0.25 ns simulation time. 

The initial velocities of carbon atoms at temperature 
T are assigned as independent Gaussian random vari- 
ables drawn from the Maxwell-Boltzmann distribution. 
All atoms are at the optimized position at i = 0. A long 
enough simulation time is used for the system to reach 
steady state. In our simulation, 5 x 10^ MD steps are 
used to ensure that the system has achieved the thermal 
equilibrium. Another 5 x 10^ MD steps are applied to 
calculate the time averaged quantities in this paper. The 
typical variation in the total energy of the system is very 
smah (< 2%). 

After we obtain the (tr^) from MD simulation, we can 
calculate the value of Young's modulus through Eq. ([7]). 
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We note that the elastic theory has been successfully ap- 
plied to describe atomic graphene system with about 400 
carbon atomst^^ In this paper, the graphene samples have 
about 200-500 carbon atoms. So we expect the Eq. ([7]) 
resulted from elastic theory can also be applicable. To 
depress the possible error created by randomness in the 
simulation, we repeat 100 independent processes for each 
value of the Young's modulus in this work. 
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FIG. 2: (Color online) The Young's modulus in graphene with 
different sizes. 

Fig. [2] shows the size dependence of the Young's mod- 
ulus. When 10 A < L < 40 A, r increases from 0.7 
TPa to 1.1 TPa with increasing size, and this value (1.1 
TPa) almost does not change with further increasing L. 
The increase of Y with increasing size also shows up in 
some studies on the Young's modulus of CNT by various 
methods, where Y increases with increasing diameter and 
reaches a saturate value i^i^i^^i^ In Fig. [2l the value of 
Y in large size sample is 1.1 TPa. This value agrees quite 
well with the recent experimental 1 ± 0.1 TPa resulti^ 
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FIG. 3: (Color online) The dependence of the Young's mod- 
ulus on temperature T for graphene with L=40 A. 

In Fig. [21 we show the temperature dependence of the 
Young's nadSflitfe in the temperature range from 100 to 



600 K. In the low temperature region [100, 500]K, Y in- 
creases for 15% as T increases. In the high temperature 
region T > 500 K, Y shows obvious decreasing behav- 
ior. This behavior indicates that the suitable tempera- 
ture region for our method is T < 500 K. If T > 500 
K, the optical phonon modes in the z direction will also 
be excited together with the flexure mode, leading to a 
larger value for the TMSVA in our MD simulation. And 
the result from Eq. ^ will underestimate the value of 
Young's modulus. 
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FIG. 4: (Color online) The isotopic disorder effect on the 
Young's modulus of graphene at T = 300 K with L = 40 A. 

Now we consider the result of the ^''C isotopic disorder 
in the pure ^^C graphene system. We expect this investi- 
gation of the isotopic disorder effect can give a useful clue 
to whether mechanical properties of graphene can be ma- 
nipulated in this way. In our simulation, to calculate the 
value of Y under a particular isotopic disorder percent- 
age, '^^C atoms are randomly substituted by certain num- 
ber of ^'^C atoms. This procedure is done independently 
in each of our 100 simulation processes for one value of 
the Young's modulus. Results are shown in In Fig.d] We 
find that the value of Y remains almost unchanged for 
the low isotopic disorder percentage (< 5%). Further in- 
creasing of the isotopic disorder percentage yields about 
15% reduction of Y. This result tells us that the purifica- 
tion of the natural graphene can not obtain a higher value 
of Y. On the other hand, about 15 % reduction of Y can 
be realized by increasing the isotopic disorder percent- 
age. However, as 20% isotopic disorder only achieves 15% 
reduction of Y, it is not an effective method to control 
the value of Y by modifying isotopic disorder percentage. 
This situation is very different from that in the thermal 
transport. The thermal conductivity has been shown to 
be very sensitive to the isotopic disorder percentage in 
the low disorder region with more than 40% reduction 
of thermal conductivity by less than 5% isotopic disor- 
der percentage; while for higher disorder percentage, the 
thermal conductivity keeps almost unchanged i^^^^^i^ So 
the thermal conductivity can be greatly enhanced by syn- 
thesizing isotopically pure nanotubesi^ 
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In conclusion, we have used MD to obtain the thermal 
vibration of graphene and then calculated the Young's 
modulus from the thermal mean-square vibration ampli- 
tude. The advantage of this approach is that we don't 
have to introduce external strain on the system, and it 
can be easily applied to study different effects on the 
Young's modulus. The theoretical results agree very well 
with the experimental ones. As an application of this 
method, we study the Young's modulus of graphene with 
different size. The temperature and isotopic disorder ef- 
fects on the Young's modulus are also investigated. It 
shows that the Young's modulus increases with increas- 
ing size when the graphene sample is smaller than 40 A, 
and reaches a saturated value in samples larger than 
40 A. The value of Y increases from 0.95 TPa to 1.1 
TPa as T increases from 100 K to 500 K. For the iso- 
topic disorder effect, Y keeps almost unchanged in the 
low disorder region (< 5%), and decreases gradually for 
15% after increasing the disorder percentage up to 20%. 
This finding provides the information that the isotopic 
disorder is not an effective method to control the Young's 
modulus of graphene. 

We should point out that why we use the constant 
value of Poisson ratio fj, = 0.17. Actually fj, also depends 
on the size, temperature and isotopic doping. By apply- 



ing external strain {sx) on the graphene in x direction, 
and using MD to record the resulted strain in y direc- 
tion (sy) under different environment, i.e. different size, 
temperatures, or isotopic doping percentage, we can ob- 
tain the value of Poisson ratio from fi — \ey/ex\- We find 
that the value of Poisson ratio will deviate from 0.17, 
which means that it will introduce some error if we use 
a constant value for Poisson ratio under all environment. 
However, as can be seen from Eq. the Poisson ratio 
appears in the expression as a factor (1 — ^^), so the er- 
ror is considerably small. For example, we find that the 
largest value for Poisson ratio is 0.22 in graphene sheet 
with L = 10 A at 300 K without isotopic doping. In 
this situation, the relative error is the largest, which is 
((l-0.222)-(l-0.172))/(l-0.172) = -2%. So through- 
out this paper, we use a constant value for Poisson ratio, 
which will introduce relative error for the Young's mod- 
ulus less than 2%. 
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